/****************************************************************************
 * Copyright (C) 2009-2010 SciTouch LLC
 * 
 * This file is part of Indigo toolkit.
 * 
 * This file may be distributed and/or modified under the terms of the
 * GNU General Public License version 3 as published by the Free Software
 * Foundation and appearing in the file LICENSE.GPL included in the
 * packaging of this file.
 * 
 * This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
 * WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
 ***************************************************************************/

#include <time.h>

#include "base_cpp/output.h"
#include "molecule/molecule.h"
#include "molecule/molfile_saver.h"
#include "graph/graph_highlighting.h"
#include "molecule/molecule_stereocenters.h"

MolfileSaver::MolfileSaver (Output &output) : v2000(false), highlighting(0),
reactionAtomMapping(0),
reactionAtomInversion(0),
reactionAtomExactChange(0),
reactionBondReactingCenter(0),
 _output(output)
{
}

void MolfileSaver::saveMolecule (const Molecule &mol)
{
   _writeHeader(_output, mol);

   if (!v2000)
   {
      _writeCtabHeader(_output);
      _writeCtab(_output, mol, false);
   }
   else
   {
      _writeCtabHeader2000(_output, mol);
      _writeCtab2000(_output, mol, false);
   }
   _output.writeStringCR("M  END");
}

void MolfileSaver::saveQueryMolecule (const Molecule &mol)
{
   if (!v2000)
   {
      _writeHeader(_output, mol);
      _writeCtabHeader(_output);
      _writeCtab(_output, mol, true);
      _output.writeStringCR("M  END");
   }
   else
   {
      if (mol.getRGroups() != 0)
      {
         time_t tm = time(NULL);
         const struct tm *lt = localtime(&tm);
         _output.printfCR("$MDL  REV  1 %02d%02d%02d%02d%02d",
            lt->tm_mon + 1, lt->tm_mday, lt->tm_year % 100, lt->tm_hour, lt->tm_min);
         _output.writeStringCR("$MOL");
         _output.writeStringCR("$HDR");
      }

      _writeHeader(_output, mol);

      if (mol.getRGroups() != 0)
      {
         _output.writeStringCR("$END HDR");
         _output.writeStringCR("$CTAB");
      }

      _writeCtabHeader2000(_output, mol);
      _writeCtab2000(_output, mol, true);

      if (mol.getRGroups() != 0)
      {
         int i, j;

         _writeRGroupIndices2000(_output, mol);

         MoleculeRGroups &rgroups = *mol.getRGroups();
         int n_rgroups = rgroups.getRGroupCount();

         for (i = 0; i < n_rgroups; i++)
         {
            const RGroup &rgroup = rgroups.getRGroup(i);

            if (rgroup.fragments.size() == 0)
               continue;

            _output.printf("M  LOG  1 %3d %3d %3d  ", i + 1, rgroup.if_then + 1, rgroup.rest_h);

            QS_DEF(Array<char>, occ);
            ArrayOutput occ_out(occ);

            _writeOccurenceRanges(occ_out, rgroup.occurence);

            for (j = 0; j < 3 - occ.size(); j++)
               _output.writeChar(' ');

            _output.write(occ.ptr(), occ.size());
            _output.writeCR();
         }

         _output.writeStringCR("M  END");
         _output.writeStringCR("$END CTAB");

         for (i = 0; i < n_rgroups; i++)
         {
            const RGroup &rgroup = rgroups.getRGroup(i);

            if (rgroup.fragments.size() == 0)
               continue;

            _output.writeStringCR("$RGP");
            _output.printfCR("%4d", i + 1);

            for (j = 0; j < rgroup.fragments.size(); j++)
            {
               const Molecule *fragment = rgroup.fragments[j];

               _output.writeStringCR("$CTAB");
               _writeCtabHeader2000(_output, *fragment);
               _writeCtab2000(_output, *fragment, true);
               _writeRGroupIndices2000(_output, *fragment);
               _writeAttachemtValues2000(_output, *fragment);

               _output.writeStringCR("M  END");
               _output.writeStringCR("$END CTAB");
            }
            _output.writeStringCR("$END RGP");
         }
         _output.writeStringCR("$END MOL");
      }
      else
         _output.writeStringCR("M  END");
   }
}

void MolfileSaver::saveCtab3000 (const Molecule &mol)
{
   _writeCtab(_output, mol, false);
}

void MolfileSaver::saveQueryCtab3000 (const Molecule &mol)
{
   _writeCtab(_output, mol, true);
}

void MolfileSaver::_writeHeader (Output &output, const Molecule &mol)
{
   time_t tm = time(NULL);
   const struct tm *lt = localtime(&tm);

   const char *dim;

   if (Molecule::haveZCoord(mol))
      dim = "3D";
   else
      dim = "2D";

   output.writeCR();
   output.printfCR("  -INDIGO-%02d%02d%02d%02d%02d%s", lt->tm_mon + 1, lt->tm_mday,
      lt->tm_year % 100, lt->tm_hour, lt->tm_min, dim);
   output.writeCR();
}

void MolfileSaver::_writeCtabHeader (Output &output)
{
   output.printfCR("%3d%3d%3d%3d%3d%3d%3d%3d%3d%3d%3d V3000",
      0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0);
}

void MolfileSaver::_writeCtabHeader2000 (Output &output, const Molecule &mol)
{
   int chiral = 0;

   if (mol.getStereocenters().haveAllAbs())
      chiral = 1;

   output.printfCR("%3d%3d%3d%3d%3d%3d%3d%3d%3d%3d%3d V2000",
      mol.vertexCount(), mol.edgeCount(), 0, 0, chiral, 0, 0, 0, 0, 0, 999);
}

void MolfileSaver::_writeAtomLabel (Output &output, int label)
{
   output.writeString(Element::toString(label));
}


void MolfileSaver::_writeMultiString (Output &output, const char *string, int len)
{
   int limit = 70;
   while (1)
   {
      output.writeString("M  V30 ");

      if (len <= limit)
      {
         output.write(string, len);
         output.writeCR();
         break;
      }

      output.write(string, limit);
      output.writeStringCR("-");
      len -= limit;
      string += limit;
   }
}


// For queries only RGroups are supported!!!!
void MolfileSaver::_writeCtab (Output &output, const Molecule &mol, bool is_query, bool is_rgroup)
{
   const Molecule *query = 0;

   if (is_query)
      query = (const Molecule *)(&mol);

   output.writeStringCR("M  V30 BEGIN CTAB");
   output.printfCR("M  V30 COUNTS %d %d 0 0 0", mol.vertexCount(), mol.edgeCount());
   output.writeStringCR("M  V30 BEGIN ATOM");

   int i;
   int iw = 1;
   QS_DEF(Array<char>, buf);
   QS_DEF(Array<int>, atom_mapping);
   QS_DEF(Array<int>, bond_mapping);

   atom_mapping.clear_resize(mol.vertexEnd());
   bond_mapping.clear_resize(mol.edgeEnd());

   for (i = mol.vertexBegin(); i < mol.vertexEnd(); i = mol.vertexNext(i), iw++)
   {
      const Atom &atom = mol.getAtom(i);
      const QueryAtom *query_atom = 0;
      int isotope = atom.isotope;
      ArrayOutput out(buf);

      if (query != 0)
         query_atom = &query->getQueryAtom(i);

      atom_mapping[i] = iw;
      out.printf("%d ", iw);

      if (atom.label == ELEM_H && isotope == 2)
      {
         out.writeChar('D');
         isotope = 0;
      }
      else if (atom.label == ELEM_H && isotope == 3)
      {
         out.writeChar('T');
         isotope = 0;
      }
      else if (query_atom == 0 || query_atom->type == 0)
      {
         if (mol.isPseudoAtom(i))
            out.writeString(mol.getPseudoAtom(i));
         else
            _writeAtomLabel(out, atom.label);
      }
      else if (query_atom->type == QUERY_ATOM_A)
         out.writeChar('A');
      else if (query_atom->type == QUERY_ATOM_Q)
         out.writeChar('Q');
      else if (query_atom->type == QUERY_ATOM_X)
         out.writeChar('X');
      else if (query_atom->type == QUERY_ATOM_RGROUP)
         out.writeString("R#");
      else
         throw Error("molfile 3000: can not save atom %d because of unsupported "
                     "query feature", i); 

      //else if (query_atom->type == QUERY_ATOM_LIST ||
      //         query_atom->type == QUERY_ATOM_NOTLIST)
      //{
      //   int k;

      //   if (query_atom->type == QUERY_ATOM_NOTLIST)
      //      out.writeString("NOT");
      //   out.writeBinaryChar('[');
      //   for (k = 0; k < query_atom->list_size - 1; k++)
      //   {
      //      _writeAtomLabel(query_atom->list[k], out);
      //      out.writeBinaryChar(',');
      //   }
      //   _writeAtomLabel(query_atom->list[query_atom->list_size - 1]);
      //   out.writeBinaryChar(']');
      //}

      int aam = 0, ecflag = 0, irflag = 0;

      if (reactionAtomMapping != 0)
         aam = reactionAtomMapping->at(i);
      if (reactionAtomInversion != 0)
         irflag = reactionAtomInversion->at(i);
      if (reactionAtomExactChange != 0)
         ecflag = reactionAtomExactChange->at(i);

      out.printf(" %f %f %f %d", atom.pos.x, atom.pos.y, atom.pos.z, aam);

      if (atom.charge != 0)
         out.printf(" CHG=%d", atom.charge);
      if (atom.radical > 0)
         out.printf(" RAD=%d", atom.radical);
      if (isotope > 0)
         out.printf(" MASS=%d", isotope);
      if (atom.explicit_valence && atom.valence != 0)
         out.printf(" VAL=%d", atom.valence);
      if (irflag > 0)
         out.printf(" INVRET=%d", irflag);
      if (ecflag > 0)
         out.printf(" EXACHG=%d", ecflag);

      if (query_atom != 0 && query_atom->type == QUERY_ATOM_RGROUP)
      {
         int k;

         MoleculeRGroups &rgroups = *query->getRGroups();

         const Array<int> &rg_list = rgroups.getSiteRGroups(i);

         if (rg_list.size() > 0)
         {
            out.printf(" RGROUPS=(%d", rg_list.size());
            for (k = 0; k < rg_list.size(); k++)
            {
               out.printf(" %d", rg_list[k] + 1);
            }
            out.writeChar(')');

            int order1 = rgroups.getAttachmentOrder(i, 0);
            int order2 = rgroups.getAttachmentOrder(i, 1);

            int n_att_points = order1 >= 0 && order2 >= 0 ? 2 : 1;

            if (n_att_points == 2 && order1 > order2)
            {
               out.printf(" ATTCHORD=(2");

               out.printf(" %d 1 %d 2", order1 + 1, order2 + 1);

               out.writeChar(')');
            }
         }
      }

      if (query_atom != 0 && query->getRGroupFragment() != 0)
      {
         MoleculeRGroupFragment &fragment = *query->getRGroupFragment();

         int val = 0;

         for (int idx = 0; idx < fragment.attachmentPointCount(); idx++)
         {
            int j;

            for (j = 0; fragment.getAttachmentPoint(idx, j) != -1; j++)
               if (fragment.getAttachmentPoint(idx, j) == i)
               {
                  val |= 1 << idx;
                  break;
               }

            if (fragment.getAttachmentPoint(idx, j) != -1)
               break;
         }

         if (val > 0)
            out.printf(" ATTCHPT=%d", val == 3 ? -1 : val);
      }

      _writeMultiString(output, buf.ptr(), buf.size());
   }

   output.writeStringCR("M  V30 END ATOM");
   output.writeStringCR("M  V30 BEGIN BOND");

   iw = 1;

   for (i = mol.edgeBegin(); i < mol.edgeEnd(); i = mol.edgeNext(i), iw++)
   {
      const Edge &edge = mol.getEdge(i);
      const Bond &bond = mol.getBond(i);
      ArrayOutput out(buf);

      bond_mapping[i] = iw;

      int bond_type = bond.type;

      if (query != 0)
      {
         const QueryBond &query_bond = mol.getQueryBond(i);
         if (query_bond.type != 0)
            bond_type = query_bond.type;
      }

      out.printf("%d %d %d %d", iw, bond_type, atom_mapping[edge.beg],
         atom_mapping[edge.end]);

      int direction = mol.getStereocenters().getBondDirection(i);

      switch (direction)
      {
         case MoleculeStereocenters::BOND_UP:     out.printf(" CFG=1"); break;
         case MoleculeStereocenters::BOND_EITHER: out.printf(" CFG=2"); break;
         case MoleculeStereocenters::BOND_DOWN:   out.printf(" CFG=3"); break;
      }

      int reacting_center = 0;

      if(reactionBondReactingCenter != 0 && reactionBondReactingCenter->at(i) != 0)
         reacting_center = reactionBondReactingCenter->at(i);

      if (reacting_center != 0)
         out.printf(" RXCTR=%d", reacting_center);

      _writeMultiString(output, buf.ptr(), buf.size());
   }

   output.writeStringCR("M  V30 END BOND");

   MoleculeStereocenters &stereocenters = mol.getStereocenters();

   if (stereocenters.begin() != stereocenters.end() || highlighting != 0)
   {
      output.writeStringCR("M  V30 BEGIN COLLECTION");

      QS_DEF(Array<int>, processed);

      processed.clear_resize(mol.vertexEnd());
      processed.zerofill();

      for (i = mol.vertexBegin(); i != mol.vertexEnd(); i = mol.vertexNext(i))
      {
         if (processed[i])
            continue;

         ArrayOutput out(buf);
         int j, type = stereocenters.getType(i);

         if (type == MoleculeStereocenters::ATOM_ABS)
            out.writeString("MDLV30/STEABS ATOMS=(");
         else if (type == MoleculeStereocenters::ATOM_OR)
            out.printf("MDLV30/STEREL%d ATOMS=(", stereocenters.getGroup(i));
         else if (type == MoleculeStereocenters::ATOM_AND)
            out.printf("MDLV30/STERAC%d ATOMS=(", stereocenters.getGroup(i));
         else
            continue;
            
         QS_DEF(Array<int>, list);

         list.clear();
         list.push(i);

         for (j = mol.vertexNext(i); j < mol.vertexEnd(); j = mol.vertexNext(j))
            if (stereocenters.sameGroup(i, j))
            {
               list.push(j);
               processed[j] = 1;
            }

         out.printf("%d", list.size());
         for (j = 0; j < list.size(); j++)
            out.printf(" %d", atom_mapping[list[j]]);
         out.writeChar(')');

         _writeMultiString(output, buf.ptr(), buf.size());
      }

      if (highlighting != 0)
      {
         if (highlighting->numEdges() > 0)
         {
            const Array<int> &h_bonds = highlighting->getEdges();
            ArrayOutput out(buf);
            out.writeString("MDLV30/HILITE BONDS=(");

            out.printf("%d", highlighting->numEdges());
            for (i = mol.edgeBegin(); i != mol.edgeEnd(); i = mol.edgeNext(i))
               if (h_bonds[i])
                  out.printf(" %d", bond_mapping[i]);
            out.writeChar(')');

            _writeMultiString(output, buf.ptr(), buf.size());
         }
         if (highlighting->numVertices() > 0)
         {
            const Array<int> &h_atoms = highlighting->getVertices();
            ArrayOutput out(buf);
            out.writeString("MDLV30/HILITE ATOMS=(");

            out.printf("%d", highlighting->numVertices());
            for (i = mol.vertexBegin(); i != mol.vertexEnd(); i = mol.vertexNext(i))
               if (h_atoms[i])
                  out.printf(" %d", atom_mapping[i]);
            out.writeChar(')');

            _writeMultiString(output, buf.ptr(), buf.size());
         }
      }
      output.writeStringCR("M  V30 END COLLECTION");
   }

   output.writeStringCR("M  V30 END CTAB");

   if (query != 0 && query->getRGroups() != 0 && !is_rgroup)
   {
      int n_rgroups = query->getRGroups()->getRGroupCount();
      for (i = 0; i < n_rgroups; i++)
         if (query->getRGroups()->getRGroup(i).fragments.size() > 0)
            _writeRGroup(output, *query, i);
   }
}

void MolfileSaver::_writeOccurenceRanges (Output &out, const Array<int> &occurences)
{
   for (int i = 0; i < occurences.size(); i++)
   {
      int occurence = occurences[i];

      if ((occurence & 0xFFFF) == 0xFFFF)
         out.printf(">%d", (occurence >> 16) - 1);
      else if ((occurence >> 16) == (occurence & 0xFFFF))
         out.printf("%d", occurence >> 16);
      else if ((occurence >> 16) == 0)
         out.printf("<%d", (occurence & 0xFFFF) + 1);
      else
         out.printf("%d-%d", occurence >> 16, occurence & 0xFFFF);

      if (i != occurences.size() - 1)
         out.printf(", ");
   }
}

void MolfileSaver::_writeRGroup (Output &output, const Molecule &query, int rg_idx)
{
   QS_DEF(Array<char>, buf);
   ArrayOutput out(buf);
   const RGroup &rgroup = query.getRGroups()->getRGroup(rg_idx);

   output.printfCR("M  V30 BEGIN RGROUP %d", rg_idx + 1);

   out.printf("RLOGIC %d %d ", rgroup.if_then + 1, rgroup.rest_h);

   _writeOccurenceRanges(out, rgroup.occurence);

   _writeMultiString(output, buf.ptr(), buf.size());

   for (int i = 0; i < rgroup.fragmentsCount(); i++)
      _writeCtab(output, *rgroup.fragments[i], true, true);

   output.writeStringCR("M  V30 END RGROUP");
}

void MolfileSaver::_writeCtab2000 (Output &output, const Molecule &mol, bool is_query)
{
   const Molecule *query = 0;

   if (is_query)
      query = (const Molecule *)(&mol);

   int i;
   QS_DEF(Array<int>, radicals);
   QS_DEF(Array<int>, charges);
   QS_DEF(Array<int>, isotopes);
   QS_DEF(Array<int>, pseudoatoms);

   radicals.clear();
   charges.clear();
   isotopes.clear();
   pseudoatoms.clear();

   for (i = mol.vertexBegin(); i < mol.vertexEnd(); i = mol.vertexNext(i))
   {
      const Atom &atom = mol.getAtom(i);
      const QueryAtom *query_atom = 0;

      if (query != 0)
         query_atom = &query->getQueryAtom(i);

      char label[3] = {' ', ' ', ' '};
      int valence = 0, isotope = 0, radical = 0, charge = 0, stereo_parity = 0,
         hydrogens_count = 0, stereo_care = 0,
         aam = 0, irflag = 0, ecflag = 0;

      if (query_atom != 0 && query_atom->type != 0)
      {
         switch (query_atom->type)
         {
            case QUERY_ATOM_RGROUP:
               label[0] = 'R';
               label[1] = '#';
               break;
            case QUERY_ATOM_A:
               label[0] = 'A';
               break;
            case QUERY_ATOM_Q:
               label[0] = 'Q';
               break;
         }
      }
      else if (atom.label == ELEM_H && atom.isotope == 2)
         label[0] = 'D';   
      else if (atom.label == ELEM_H && atom.isotope == 3)
         label[0] = 'T';
      else if (mol.isPseudoAtom(i))
      {
         label[0] = 'A';
         pseudoatoms.push(i);
      }
      else
      {
         const char *str = Element::toString(atom.label);

         label[0] = str[0];
         if (str[1] != 0)
            label[1] = str[1];

         if (atom.isotope >= -18 && atom.isotope <= 12)
            isotope = atom.isotope;

         if (isotope != 0)
            isotopes.push(i);
      }

      if (reactionAtomMapping != 0)
         aam = reactionAtomMapping->at(i);
      if (reactionAtomInversion != 0)
         irflag = reactionAtomInversion->at(i);
      if (reactionAtomExactChange != 0)
         ecflag = reactionAtomExactChange->at(i);

      if (atom.explicit_valence && atom.valence > 0 && atom.valence < 14)
         valence = atom.valence;

      if (atom.charge >= -15 && atom.charge <= 15)
         charge = atom.charge;

      if (charge != 0)
         charges.push(i);

      if (atom.radical >= 0 && atom.radical <= 3)
         radical = atom.radical;

      if (radical != 0)
         radicals.push(i);

      output.printfCR("%10.4f%10.4f%10.4f %c%c%c%2d"
                    "%3d%3d%3d%3d%3d"
                    "%3d%3d%3d%3d%3d%3d",
         atom.pos.x, atom.pos.y, atom.pos.z, label[0], label[1], label[2], 0,
         0, stereo_parity, hydrogens_count, stereo_care, valence,
         0, 0, 0,
         aam, irflag, ecflag);
   }

   for (i = mol.edgeBegin(); i < mol.edgeEnd(); i = mol.edgeNext(i))
   {
      const Edge &edge = mol.getEdge(i);
      const Bond &bond = mol.getBond(i);

      int stereo = 0;
      int topology = 0;
      int reacting_center = 0;

      int direction = mol.getStereocenters().getBondDirection(i);

      switch (direction)
      {
         case MoleculeStereocenters::BOND_UP: stereo = 1; break;
         case MoleculeStereocenters::BOND_DOWN: stereo = 6; break;
         case MoleculeStereocenters::BOND_EITHER: stereo = 4; break;
      }

      int bond_type = bond.type;

      if (query != 0)
      {
         const QueryBond &query_bond = query->getQueryBond(i);
         if (query_bond.type != 0)
            bond_type = query_bond.type;
      }

      if(reactionBondReactingCenter != 0 && reactionBondReactingCenter->at(i) != 0)
         reacting_center = reactionBondReactingCenter->at(i);

      output.printfCR("%3d%3d%3d%3d%3d%3d%3d", edge.beg + 1, edge.end + 1, bond_type,
         stereo, 0, topology, reacting_center);
   }

   if (charges.size() > 0)
   {
      output.printf("M  CHG%3d", charges.size());
      for (i = 0; i < charges.size(); i++)
         output.printf(" %3d %3d", charges[i] + 1, mol.getAtom(charges[i]).charge);
      output.writeCR();
   }

   if (radicals.size() > 0)
   {
      output.printf("M  RAD%3d", radicals.size());
      for (i = 0; i < radicals.size(); i++)
         output.printf(" %3d %3d", radicals[i] + 1, mol.getAtom(radicals[i]).radical);
      output.writeCR();
   }

   if (isotopes.size() > 0)
   {
      output.printf("M  ISO%3d", isotopes.size());
      for (i = 0; i < isotopes.size(); i++)
         output.printf(" %3d %3d", isotopes[i] + 1, mol.getAtom(isotopes[i]).isotope);
      output.writeCR();
   }

   for (i = 0; i < pseudoatoms.size(); i++)
   {
      output.printfCR("A  %3d", pseudoatoms[i] + 1);
      output.writeString(mol.getPseudoAtom(pseudoatoms[i]));
      output.writeCR();
   }
}

void MolfileSaver::_writeRGroupIndices2000 (Output &output, const Molecule &query)
{
   int i, j;

   QS_DEF(Array<int>, atom_ids);
   QS_DEF(Array<int>, rg_ids);

   atom_ids.clear();
   rg_ids.clear();

   for (i = query.vertexBegin(); i < query.vertexEnd(); i = query.vertexNext(i))
   {
      if (query.getQueryAtom(i).type != QUERY_ATOM_RGROUP)
         continue;

      const Array<int> &rg_list = query.getRGroups()->getSiteRGroups(i);

      for (j = 0; j < rg_list.size(); j++)
      {
         atom_ids.push(i + 1);
         rg_ids.push(rg_list[j] + 1);
      }
   }

   if (atom_ids.size() > 0)
   {
      output.printf("M  RGP%3d", atom_ids.size());
      for (i = 0; i < atom_ids.size(); i++)
         output.printf(" %3d %3d", atom_ids[i], rg_ids[i]);
      output.writeCR();
   }
}

void MolfileSaver::_printOccurenceRanges (Output &out, const Array<int> &occurences)
{
   for (int i = 0; i < occurences.size(); i++)
   {
      int occurence = occurences[i];

      if ((occurence & 0xFFFF) == 0xFFFF)
         out.printf(">%d", (occurence >> 16) - 1);
      else if ((occurence >> 16) == (occurence & 0xFFFF))
         out.printf("%d", occurence >> 16);
      else if ((occurence >> 16) == 0)
         out.printf("<%d", (occurence & 0xFFFF) + 1);
      else
         out.printf("%d-%d", occurence >> 16, occurence & 0xFFFF);

      if (i != occurences.size() - 1)
         out.printf(", ");
   }
}

void MolfileSaver::_writeAttachemtValues2000 (Output &output, const Molecule &fragment)
{
   if (fragment.getRGroupFragment() == 0)
      return;

   MoleculeRGroupFragment &rg_fragment = *fragment.getRGroupFragment();

   if (rg_fragment.attachmentPointCount() == 0)
      return;

   RedBlackMap<int, int> orders;
   int i;

   for (i = 0; i < rg_fragment.attachmentPointCount(); i++)
   {
      int j = 0;
      int idx;

      while ((idx = rg_fragment.getAttachmentPoint(i, j++)) != -1)
      {
         int *val;

         if ((val = orders.at2(idx + 1)) == 0)
            orders.insert(idx + 1, 1 << i);            
         else
            *val |= 1 << i;
      }

   }

   output.printf("M  APO%3d", orders.size());

   for (i = orders.begin(); i < orders.end(); i = orders.next(i))
      output.printf(" %3d %3d", orders.key(i), orders.value(i));

   output.writeCR();
}
